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simulated.  Turbulence  Is  introduced  using  a two  layer  mixing  length  model 
appropriate  to  curved  three-dimensional  wall  jets.  Typical  results  quantifying 
Jet  spreading,  Jet  growth,  nominal  separation  and  Jet  shrink  effects  due  to 
cross  flow  are  presented, 

A class  of  cases  was  l^estigated  roughly  possessing  initial  flow 
angularity  and  adverse  pressure  gradients  prototypic  of  those  on  the  blown 
surfaces  of  typical  propulsive  lift  systems  such  as  the  havy/Rockwell  XFV-12A 
thrust  augmented  wing.  Results  obtained  from  the  computational  model  Indleata 
that  if  the  initial  total  velocity  is  kept  fixed,  then  the  Introduction  of  the 
cross  flow  enhances  the  freestream  decay  rate  of  the  peak  of  the  velocity 
component  in  the  freestream  direction.  In  addition,  the  entrainment  quantity 
and  its  rate  decrease  with  Increased  cross  flow.  The  three-dimensional 
phenomena  not  only  influence  the  effect  of  taper  on  the  boundary  layer  control 
characteristics  of  a Coanda  flap,  but  also  indicate  a "Jot  shrink"  which  could 
be  a mechanism  promoting  end-wall  separation.  To  our  knowledge,  our  model  is 
the  first  to  quantify  such  trends.  Both  should  be  considered  In  the  design  of 
any  propulsive  lift  system.  Finally,  the  effect  on  the  prescribed  external 
adverse  pressure  gradient  in  the  presence  and  absence  of  cross. flow  has  also 
been  examined.  From  the  limited  results,  the  spsnwise  separation  line  moves 
progressively  further  upstream  with  Increasing  cross  flow. 
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FOREWORD 


This  docunant  dflicrlbaa  computational  atudlaa  of  thraa-dlmanalonal 
Incompraaalbla  turbulant  wall  jata  with  large  croas  flow.  Tha  effort 
waa  parformad  during  tha  period  September  19,  1977  to  Saptaobar  18.  1978 
and  waa  aponaored  by  tha  Naval  Air  Development  Canter  under  Contract 
N62269-77-C-0412. 

Tha  technical  monitor  for  thla  atudy  waa  Dr.  K.A.  Graen. 
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ABSTRACT 


The  flow  field  of  three  dimensional  Incompressible  wall  jets 
prototypic  of  thrust  augmenting  ejectors  with  large  cross  flow  Is  solved 
using  a very  efficient  centared-Euler  scheme  in  an  orthogonal  curvilinear 
coordinate  system.  The  computational  model  treats  Initial  conditions 
with  arbitrary  velocity  profllee  at  or  downstream  of  the  Jet  exit. 

An  averaging  approach  la  employed  for  the  first  few  marching  steps  to 
overcome  spurious  numerical  oaclllatlona  associated  with  arbitrary 
initial  profiles.  Laminar  aa  well  aa  turbulent  wall  jats  are  simulatad. 
Turbulence  la  introduced  uaing  a two  layer  mixing  length  modal  appropriate 
to  curved  thrae-dimanaional  wall  jats.  Typical  rasulta  quantifying  Jet 
apreading,  jet  growth,  nominal  aeparatlon  and  Jet  shrink  effects  due  to 
cross  flow  are  presented. 

A class  of  ceaas  was  Invsatlgatsd  roughly  poesaealng  initial  flow 
angularity  and  adverae  preaaura  gradlanta  prototyplc  of  Choae  on  the 
blown  surfaces  of  typical  propulsive  lift  systems  such  as  the  Navy/ 

Rockwell  XFV-12A  thrust  augmented  wing.  Results  obtained  from  the 
computational  model  Indicate  that  If  the  Initial  total  valoclty  Is  kept 
fixed,  then  the  Introduction  of  the  cross  flow  enhancas  tha  fraastraam 
decay  rate  of  the  peak  of  tha  valoclty  componant  in  the  freestreem 
direction.  In  addition,  the  entrainment  quantity  and  Its  rate  daersaee 
with  increased  cross  flow.  The  three-dlmenelonal  phenomena  not  only 
influence  the  effect  of  taper  on  the  boundary  layer  control  cheractaristlcs 
of  a Coanda  flap,  but  alao  indicate  a "jet  shrink"  which  could  be  a 
mechanism  promoting  and-wall  separation.  To  our  knowledge,  our  model 
la  tha  flrat  to  quantify  such  trends.  Both  should  be  considered  in  the 
design  of  any  propulelve  lift  system.  Finally,  the  effect  on  the 
prescribed  externsl  adverse  pressure  gradient  In  tha  presence  and  absence 
of  cross  flow  lus  also  been  examined,  From  the  limited  results,  the 
spanwlss  separation  line  movaa  progressively  further  upatreem  with 
Increasing  cross  flow. 
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NOMENCLATURE 


b 

f.g 

Ki.Kj 


M.N.Q.R 


Q 


u,v,w 


u'v',  v'w' 

Si 

S 

(*iy**) 


y* 

h 


4* 

e 

<1,K2 

'll 


Defined  In  Eq.  (14a) 

Reduced  vector  potentials  (Eqa.  (8a)  and  (8b)) 

Metric  coefficients 
Geodesic  curvatures 
Constants  defined  on  page  1-12 
Constanta  defined  on  page  1-12 
Entrainment  quantity 

Xiy,2  components  of  velocity  (see  Fig.  3) 

Reynolds  stress  correlations 
Arc  length  along  x coordinate 
Area  of  Integration 

Orthogonal  curvilinear  coordinates  parallel  and 
perpendicular  to  wall  (see  Fig.  3) 

Height  of  first  turbulent  layer  in  two  layer  eddy 
viscosity  model 

Constant  appearing  in  eddy  viscosity  model,  page  1-9 

Eddy  viscosity  (page  1-9) 

Reduced  eddy  viscosity  (page  1-9) 

Reduced  addy  viscosity  defined  In  Eq.  (14a) 

Radius  of  curvature  of  z-constant  and  x«conatant  linns 
on  wall  Jet  surface 

Component  of  vector  potential  (page  1-8) , velocity 
potential  (page  1-20) 

Component  of  vector  potential  (page  1-8) 
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NOMENCLATURE  (Continued) 


0 Source  strength  (page  1-20) 

V Kinematic  viscosity 

T Shear  stress  ■ p(e+v)  ^0*  + 

d Initial  jet  flow  angle  parameter,  Fig.  8 

Dummy  variables,  page  1-20 

Subscripts 

e External  flow 

l,j  Modal  locations  for  difference  operators 

(Mote  also  that  subscript  notation  is  used  in  what 
follows  to  denote  partial  differentiation.) 

Superscript 


n 


Time  level  and  nodal  location  for  difference  operators 


PART  I - THEORETICAL  ANALYSIS  AND  RESULTS 


1.0  Introduction 

Modem  naval  aircraft  can  reduce  strike  force  vulnerability  through 
greater  dlapersal.  A way  cf  achieving  this  dispersal  would  be  the 
davelopmeut  and  deployment  of  vertical  lift-off  and  landing  aircraft 
operatiug  from  ships  that  are  significantly  smaller  than  the  currently 
operating  carriers.  One  propulsion  concept  designed  to  achieve  a 
vertical  and  short  take-off  and  landing  (V/STOL)  aircraft  Involves  the 
use  of  thrust  augmenting  ejectors  to  amplify  the  thrust  of  the  engines 
in  the  vertical  mode.  This  tech.iiique  is  utilized  in  the  current  XFV-12A 
aircraft,  To  achieve  adequate  accelerations  for  typical  payloads,  a 
high  augmentation  ratio  ((^)  is  required.  Various  augraenter  designs  have 
been  proposed  to  achieve  these  high  d)  values  as  well  as  integrate  well 
with  the  aircraft  structure.  In  the  XFV-12A,  an  ejector  system  composed 
of  a centerbody  and  two  Coanda  wall  jets  is  currently  under  development. 

A central  feature  of  the  flow  fields  produced  by  this  device  is  three 
dimensionality.  This  has  been  particularly  evident  in  subacale  flow 
visualization  on  the  Coanda  surfaces.  It  is  believed  chat  these  flow 
processes  nay  be  important  toward  d maximization.  The  augmentation  ratio 
depends  on  the  n»menCun  flux  and  the  dynamic  head  losses  downstream  of 
the  Coandas  and  the  centerbody  nozzle,  and  in  particular,  in  the  diffusor 
section  of  the  augmenter.  Maintenance  of  a high  momentum  and  minimiza- 
tion of  extraneous  motions,  while  accelerating  as  quickly  as  possible 
the  secondary  flows  using  the  primary  Jets,  is  of  paramount  importance. 

An  example  of  an  extraneous  motion  is  associated  with  certain  vortex  struc- 
tures which  have  been  observed  in  augmenter  wing  configurations.  It  has  been 
postulated  that  these  motions  and  general  flow  turning  are  induced  by 
the  three-dimensional  features  of  the  configuration.  If  the  wing  flap 
nn  which  the  Coandas  are  mounted  has  a tapered  diffusor,  swept  back 
trailing  edge,  and  small  aspect  ratio  planform,  considerable  cross  flow 
in  the  secondary  can  arise.  Under  these  circumstances,  this  turning  is 
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further  enhanced  by  spanwlse  preasure  gradients  arising  along  the  tapered 
canterbody  and  Coanda  slots.  These  factors  can  also  lead  to  considerable  flow 
nonunlfomltles  both  in  the  apanwlse  and  chordwise  directions.  To 
compound  the  situation,  flow  surveys  for  three  dimensional  research 
confl^.uratlons  of  this  type  have  Indicated  regions  of  separation  near 
the  andwslla.  The  amount  and  nature  of  blowing  In  these  regions  to 
maintain  attached  flow  and  craate  the  highest  posslbls  momsntum  flux 
across  some  downatraam  control  surface  needs  quantification. 

From  the  foregoing  cone iderat Ions,  it  la  obvloua  that  an  analysis 
which  asBumas  mean  flow  quantities  such  as  that  employed  in  tha  usual 
one-dimenaional  derivations  can  be  considerably  In  error.  In  fact,  if 
it  is  aasumad  in  such  a treatment  that  tha  flow  is  turned  by  an  angle 
u uniformly  along  the  apan,  the  reduction  in  <t>  Is  only  cos  u),  which  is 
considerably  less  than  rough  comparisons  between  quast-two-dlmenalonal  con- 
figurations and  actual  three-dimensional  arrangements  would  Indicate.  One 
way  of  understanding  the  foregoing  relationships  is  through  thsoratical 
modeling  which  can  provide  a means  of  reducing  the  high  cost  of  powered 
lift  testing.  Unfortunately,  existing  methodology  has  been  llraltad  in 
the  past  to  two-dimensional  flows  for  the  analysis  of  wall  Jets  and 
complete  ejector  systems.  Analytical  methods  and  computational  algorithms 
are  therefore  necessary  to  compute  three-dimensional  flows  typical  of 
reality. 

To  shed  light  on  typical  flow  patterns  encountered  due  to  the  effect 
of  taper  and  sweep  on  augmenter  wings  as  well  as  upper-surface-blown 
configurations,  a study,  "Three-Dimensional  Flow  of  a Wall  Jet,"  was 
initiated  by  the  Naval  Air  Development  Canter  to  investigate  wall  Jet  flows 
which  exemplify  typical  features  of  more  complex  propulsive  lift  applications. 
The  purpose  of  this  study  has  been  to  apply  modern  computational  methods  to 
the  treatment  of  wall  Jet  flows  with  three  dimensionality.  In  a previous 
phase  of  the  effort,  small  cross  flow  wall  Jets  were  considered.*  This 
report  relaxes  that  assumption  and  considers  large  cross  flows  that  occur 
in  practice. 
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The  formulation  employs  boundary  layer  equations  in  an  orthogonal 
curvilinear  coordinate  system.  If  the  distance  from  the  jet  exit  is 
sufficiently  large  to  establish  complete  mixing,  the  Jet  exit  height  is 
small  compared  to  a characteristic  radius  of  curvature,  and  the  Reynolds 
number  based  on  the  exit  height  Is  large,  the  ordsr  of  magnitude  analysis 
given  in  Ref.  1 extended  to  three  dimensions  Indicates  that  the  wall  jet 
equations  are  substantially  the  same  as  the  three-dimensional  boundary 
layer  aquations.  The  basic  idea  In  both  the  wall  Jet  and  boundary  layer 
approximations  is  the  same,  namely,  that  diffusional  gradients  for  the 
vortlclty  in  the  direction  normal  to  the  surface  are  dominant  terms  in  the 
equations.  In  fact,  half  of  the  wall  Jet  has  a boundary  layer  character 
due  to  the  no-sllp  condition  on  the  surface.  Furthermore,  to  dominant 
order,  the  pressure  is  independent  of  the  coordinate  normal  to  the  surface. 
Equilibration  of  centrifugal  force  between  Che  streamlines  with  pressure 
gradients  across  them  is  accomplished  with  Che  second-order  pressure  term 
and  is  expressed  in  the  second-order  momentum  aquation  in  the  normal 
direction.  The  complete  mixing  condition  is  synonymous  with  merger  of 
the  boundary  layer  from  the  wall  and  the  shear  layer  from  the  Jet  exit. 
Assuming  that  the  characterlscic  distance  for  merger  D is  determined  by  a 
spreading  angle  X of  the  order  of  2.86*  appropriate  for  submerged  turbulent 
free  shear  layers  measured  by  Reichardt,^  then  the  merger  distance  la 
d ctn  X where  d is  the  slot  height.  For  an  aircraft  similar  to  the  XFV-12A, 
assuming  that  d 'V'  3 cm  and  the  wing  chord  L 300  cm,  this  estimate  shows 
that  D/L  0,2.  Note,  however,  that  this  can  be  reduced  to  practically 
zero  if  the  wall  boundary  layer  of  the  flow  emanating  upstream  of  the 
nozzle  exit  completely  fills  the  exit  location. 

For  the  computational  model,  a transformation  is  incorporated  to 
stretch  the  coordinate  normal  to  the  flow.  At  streamwisa  planes,  the 
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resulting  nonlinear  partial  differential  equations  are  treated  aa 
ordinary  differential  equations,  incorporating  source  terms  Involving 
partial  derivatives  representing  the  upstream  history  of  the  flow  field. 

Tliese  are  solved  using  a very  efficient  two-point  boundary  value  finite- 
difference  method  devised  by  Keller  and  Cebacl’’'’  known  as  the  "box  method," 
In  the  code,  the  turbulence  Is  Introduced  using  a two  layer  mixing  length 
model  appropriate  to  three-dimensional  wall  Jets. 

2.0  Problem  Deecriptlon  and  Formulation 

To  fix  the  Ideas  regarding  three-dimensional  wall  jet  flows  i event 
to  propulsive  lift  applications,  consider  two  prototyplc  configurations  as 
Indicated  In  Figures  1 and  2.  The  generic  arrangement  shown  in  Figure  1 Is 
relevant  to  an  augmenter  ejector  wing  of  the  type  used  on  the  XFV-12A 
without  the  csnterbody  nozzle.  In  Figure  2,  the  shape  Indicated  corresponds 
to  an  application  Involving  boundary  layer  control  and  superclrculation 
development.  For  both  cases,  the  development  of  the  wall  jet  over  the 
curved  surfaces  Is  of  Interest.  In  the  augmenter  of  Figure  1,  Coande 
jets  flow  over  the  surfaces  and  $2  emanating  from  slots  and  T2  end 
provide  a vertical  lift  force.  In  the  analysis  of  the  flow  field  over  S^, 
we  suppress  the  Influence  of  the  surfaces  S2,  A,  and  C.  In  addition,  we 
assume  that  the  secondary  flow  produced  by  entrainment  resulting  from  the 
primary  Jets  emanating  from  T^^  and  T2  is  known  3:  pjrtori.  In  actuality, 
these  must  be  computed  as  an  Integral  part  of  the  problem.  For  tractabillty, 
we  restrict  our  consideration  to  the  Indicated  (parabolic)  formulation 
since  it  la  a building  block  to  a later  analysis  of  the  primary  secondary 
Interaction.  For  the  configuration  of  Figure  2,  the  orientation  Is  similar, 
and  we  neglect  the  primary  secondary  Interaction  features.  Thus,  both 
configurations  lead  to  the  problem  of  the  development  of  a 3-D  wall  Jet 
over  a curved  surface  in  which  boundary  conditions  are  apecifled  on 
some  Interfacial  layer  with  the  external  or  entrained  flow.  It  should 
be  noted  that  the  mixing  downstream  of  the  slots  is  with  a flow  above 
the  slot  which  obeys  the  no  slip  condition  at  the  slot  trailing  edge. 
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In  Figure  3,  the  prototyplc  geometry  for  the  flow  over  the  aurfacea 
Is  ahown  with  the  prescribed  Initial  and  boundary  velocity  distributions 
schematically  Indicated,  k curvilinear  coordinate  system  Is  used  which 
consists  of  conveniently  oriented  lines  In  the  surface  and  normals  to 
them  as  shown  In  the  figure.  The  Initial  distribution  of  velocities  la 
specified  along  some  line  x^XQ^constant  say,  which  tnay  be  the  Jet  exit. 

In  the  wall  Jet/boundary  layer  approximation  the  appropriate  aquations 
describing  the  wall  Jet  flow  arsi 

Continuity 

0/3x)^h2u)  + 13/3*) (h^w)  + (3/3y) (h^h^v)  - 0 . U) 


x-Monentum 


u 3u  , w 3u  . 3u 
hj^  3x  h2  3*  ^ 3y 


uwK^  + w^Kj 


s-Momentum 


_1 

ph 


3y 


u 3v 
hj^  3x 


w 3w 
h.^  3* 


V 


3y 


“ uwK^ 


+ u*K^ 


1 3p  ( 3 
ph.^  3y 


3w 

3y 


Hare,  h^  and  hj  are  metric  coefficients  and  are  functions  of  x and  z,  and 
the  parametsrs  and  are  known  as  the  geodesic  curvatures  of  the  curves 
z ■■  constant  and  x constant,  respectively. 

The  boundary  conditions  for  Eqs.  (.1)  through  (3)  for  zero  mass  transfer 
at  the  wall  and  compatibility  with  the  external  flow  are 


y ■ 0 


u,w,v  ■ 0 


y uo 


u u (x,z) 


1-6 


w -*■  w„lx,z)  . (ft) 


NADC-77163-30 


At  Indicated  earlier,  the  previous  equations  are  transformed  by 
defining 


X ■ X z 


1/2. 


- r n ■ (u^/vs^)  y 


(5) 


and  introducing  a two-component  vector  potential  given  by 


hu.ii  hw-^  hhv--(^  + ^\. 
"2“  3y  *'1''  3y  "12'^  \3x  3*/  ‘ 


(6) 


Here  e^,  which  denotes  the  ere  length  along  the  x coordinate,  is  defined  by 


'1  ” r "i’*’'  • 


(7) 


In  addition,  the  dltnenalonless  variables  f end  g related  to  ij)  and  1))  are 
defined  by 


l/2v 


tp  ■ (u^vij^)  hjfCx,*.*!) 


(Be) 


,1/2 


- (u^VSj^)  hj^(w^/u^)g(x,2,n)  • 


(8b) 


The  parameters  and  in  Eqs.  (2)  and  (3)  are  defined  by 


1 ‘‘^l 

K ■<  - — = 

1 h^hj  dz 


end 


^ 

^2  h^hj  dx  • 


3.0  Eddy  Viscosity  Modal 

Equations  (2)  and  (3)  contain  the  Reynolds  shear  stress  terms  -u' v' 


and  -v^V*^ . These  quantities  must  be  further  characterized  in  order  to 
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solva  Eqs.  (l)-(3).  To  achl«v«  this  closure,  suitable  eddy  viscosity  models 
are  employed.  In  two  dimensions,  such  simulations  have  been  successfully  used 
by  Dvorak,"  Ramaprlsn,’  and  Wilson  and  Goldstein.®  On  a heuristic  basis  subject 
to  future  experimental  validation,  extension  to  three-dimensional  flows  was 
effactad  analogous  to  an  approach  used  for  the  same  purpose  in  connection  with 
boundary  layer  flows  in  Ref.  5.  For  the  well  jet  case,  we  assume  that  the  curve' 
tura  effects  era  modeled  in  terns  of  the  principal  curvatures  in  the  x and  t 
directions.  With  this  viewpoint,  the  Reynolds  shear  stress  tarns  are  assumed 
to  be  In  the  form 


•u'v' 


cu  ■ G 

y 


uK. 


1 - 


2 1 


l'fK2y  u, 


yj 


-V  W ■ Gw 


wK 


1 - 


1_ 

l+KlX  Wy 


^ • 


(9) 


The  quantity  e in  Eq.  (9)  is  the  eddy  viscosity  of  the  corresponding 
plane  flow.  It  la  assumed  to  be  tns  same  In  both  the  x and  t dirsetione 
and  is  represented  by  a two-layer  model.  The  second  term  inside  the  bracket 
in  Eq.  (9)  Is  due  to  curvature  where  and  Kj  denote  the  radius  of  curvature 
of  i»conatant  and  x<>constant  lines.  Similar  curvature  terms  for  two- 
dimensional  wall  jets  have  been  used  by  Dvorak*  and  Wilson.'  Referring 
to  Figure  4,  the  structure  of  the  two-laysr  eddy  viscosity  modal  is  sa 
follows: 


First  Layer 

F • <0.435  y)*  /u*  + w*  0 < y < y* 

Second  Lavar 

G ■ <0.125  Vuy  + w*  y y* 

where  at  y ■ yj^ 


1 >/  U*  + wj  - >/u'  + w^  1 


>■  <10) 
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and  y*  is  obtained  by  Imposing  continuity  In  e at  y**y*.  This  yields 
y*  ■ g'—jg  yj^.  Two-layer  eddy  viscosity  models  similar  to  Eq.  (10)  have 
been  employed  by  Ramaprian^  for  two-dimensional  wall  jets  which  give 
good  comparisons  with  experiments.  In  these  cases,  the  term  appearing 
Inside  the  square  root  In  Eq.  (10)  le  zero. 

With  the  concept  of  eddy  viscosity  and  with  the  previous  trensformed 
variables  I It  can  be  shown  that  the  system  of  Eqs.  (1)  through  (4)  can 
be  written  as 


x-Moaentum  Equation 

(bf)'  + P.ff"  -f  P,[l  - (f')M  + P.[l  - f#'] 


6 

xP 


+ p.f'g  + Poll-  - (g’)M 


10 


z "Momentum  Equation 


ill 

« fit  M 4.  p /g 

. ii  . 

f"  ii.\ 

ax 

ax  7\* 

dz 

3z  / 

(bg")'  + P.fg"  + P,  (1  - £*g')  + PJl  - (g')*] 


+ P^gg"  + Pg(l  - (f')*l 


- xP 


MSI  W'4««WM>«e,irr 
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Hare  the  primes  denote  differentiation  with  respect  to  rii  end 

b “ l+e"^  ■ e/v  f - u/u^  g'  w/w^  . (I4a) 

The  coefficients  to  Pj^q  are  functions  of  u^,  w^,  h2«  K^»  and  K2 

and  are  given  by  the  following  formulas t 


P^  - (M+l)/2  - Sj^Kj 


P2-M 


‘1*^2 


^ ■ -A> 


M 


u h,  3x 

e 1 


N 


u h,  3z 

0 L 


(Ub) 


8,  3w  s,  3w 

0 , L.  — • B « 1 _S. 

^ u^hj^  3x  u^h2  3a  ' 

In  order  to  solve  eqs.  (11)  through  (13),  Initial  conditions  are 
required  at  a starting  plane,  In  the  case  of  the  boundary  layer  problem, 
the  Initial  conditions  at  x*0  and  z*0  planes  are  obtained  by  solving  the 
limiting  form  of  Eqs,  (11)  and  (12).  For  a wall  Jet,  Initial  velocity 
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profiles  are  praacrlbed  at  some  dovmstream  plane  and  along  the  z-0  plane, 

"attachment  line"  equations  are  solved.  The  attachment  line  equations  are 
obtained  by  differentiating  the  z-momentum  equation  with  respect  to  z and 
setting 


3c  3*  3i  3,»  ^ • 


The  resulting  attachment  line  equations  valid  at  the  2*0  plane  are 

(bf")'+  P^ff"  + P^tl  " (f)*]  + Pjgf"  - xPio[''  ^ 1^] 


(bg")'  + P^fg"  + P^(l  - f'g')  + Pjtl  - (g')»] 


P as"  ■ xP  [f ' - s" 

^3**  *^10  [ 3x  * 3xJ  * 


Here,  g'  la  defined  ae  w./w  , and  its  definition  corraaponde  to  L'Hoapltal's 

S 9 

s 

rule  applied  to  the  expression  for  g*  used  previously. 

4 . 0 Finite  Difference  Equations 

First,  reduce  the  system  (11)-(12)  to  the  first  order  system 


f - u 


u'  • V 


g'  ■ w 


w'  - t 


1,.;  .. 
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(bv)'  + P^fv  + P^Cl  - u*)  + Pj[l  - uw] 


+ PfiVg  + P„[l  - w*]  - “ Is  - + P7  R If)]  ^21; 


(bt)'  + Pj^ft  + P^(l  - uw)  + Pj(l  - w*)  + Pj^gt 


+ Pg(l  - u^)  . xP^Q  sS  - t 


ii  + p L Iw  . t i&'ll  . 

3x  3z  Sz/J 


Using  th«  notation  of  Ref.  1 aasocloitad  with  the  box  method  described  there 


In  Refs.  3,  4,  5, 

and 

9,  we 

let 

X.  ••  constant 

X 

■ * , 

+ k 

n ■ 

1,2,. 

..,N 

0 

n 

n-1 

n 

Zq-O 

*1 

“ *1-1 

1 - 

1,2 , • 

...I 

0 

1 

0 

• 

ZJ 

I 

J - 

CM 

s-^ 

1 • 1 <J 

Than,  using  the  box  method,  we  have 


.n,i  -n,l 
!l— ■ ^n.l 
h^  ''j-1/2 


n , 1 n , 1 

II -..“j-l  ■ v".‘ 

hj  "J-l/J 


n,l  n,l 

fj liL-i.  ■ 

h^  ''j-1/2 


n,l  n,l 
w/  -•  w/. 


;i-i/2  • 


. 'A-V-l- 
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We  use  the  notation 


--  _ n-1/2, 1-1/2  ^ n-1/2 

J-1/2  ■ i-i/2 


and 


— 1 / n,l  , n,l-l  , 

^ - 2 ^ '^J-1/2  ] 

— 1 / n,l  , n-l,l\ 

\ - 2 ('^j-1/2  + '^j-1/2)  • 

Equation  (21)  becomes,  with  the  box  centered  at 

- - * **2^^  ■ '^j-1/2^  ■ ^5^^  ” “j-l/2''j-l/2^ 

" ^6^''®^j-l/2  ■ ^8*'^  “ ''j-1/2^ 


Vl/2^10 


n n-i 


'j-1/2 


"j-1/2 


n n-1 


+ P-, 


w 


J-1/2 


- ^1-1^  - 

r,  " ''j-1/2  r 


^®1  “ *1-1^ 


1 


(27) 


Equation  (22)  and  the  attachment  line  equations  (13)-(1^)  are  discretized 
similarly.  Details  of  the  procedure  are  given  in  Ref.  9. 
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The  solution  procedures  Involve  the  following  steps: 

(1)  Solve  the  attachment  line  equations  (15)-(16)  with  boundary 
conditions  (13)  at  x “ and  z - 0 assuming  Initial 
conditions  on  x ■ Xq* 

(2)  March  in  the  z-direction  along  the  plane  x ■ and  solve 
equations  (17)- (22)  with  boundary  conditions  (13)  for  the 
unknowns  (f ,u,v,g«w,t) . 

(3)  Repeat  steps  (1)  and  (2)  for  the  next  x-plane,  x ■ X21  and 
so  on. 


The  most  efficient  way  to  solve  the  finite  difference  equations  is 
to  use  a pseudo-Newton's  relaxation  scheme.  These  equations  may  be  written 
as  a system  of  nonlinear  algebraic  equations  by  writing 

<t(u)  ■ 0 


where 


,-n,i  n,i  n,i  n,l  n,i 

J ■ (fj  . , Vj-  , 8j  , . tj  . 


Then,  the  relaxed  Newton's  method  becomes 


1 

I 


(28a) 
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The  method  Is  said  to  have  converged  when 


liSu^^  ^^11  < e (a  prescribed  error  tolerance) 


We  call  Eq.  (28)  a pseudo-Newton's  method  because  we  linearize  the 
b terms  in  Eqa.  (21)  and  (22)  by  evaluating  their  at  v-2  before  computing 
the  Jacobian  matrix,  3<t>/3u.  Consequently,  this  algorithm  will  not  be 
quite  quadratlcally  convergent.  We,  therefore,  employ  relaxation  (oi  ^ 1) 
to  accelerate  it.  Remarkably,  underrelaxation  (u  < 1)  works  very  well, 
while  overralaxation  (ui  > 1)  diverges.  Values  of  c)  of  O.S,  0.6,  0.7, 

0.3,  and  0,9  all  give  good  results  with  u •*  0.7  found  to  be  the  overall 
beat  for  some  of  our  computational  experiments . 

An  Important  feature  of  Keller's  box  method  la  that  the  Jacobian 
matrix  can  be  put  Into  block  trldlagonal  form  and  vary  efficient  elimi- 
nation schemes  can  be  employed  for  solving  Eq.  (28a). 

Minor  Difficulties  with  the  Kumerical  Algorithm 

When  starting  at  x ■ 1 with  supplied  completely  merged  wall  jet 

velocity  profiles  described  In  Fig.  5,  unnatural  oscillations  developed 

In  the  solution.  This  difficulty  was  eliminated  completely  by  employing 

the  following  "trick."  The  first  10  mesh  points  in  the  x-directlon  were 
-4 

set  at  k ■ 10  . For  the  first  five  planes  In  the  x-direction  and  all 

n 

points  in  Che  z-directlcn  in  these  planes,  an  average  value  was  used  for 
past  points,  i.e., 


CAN  BE  PRESCRIBED  AS  A FUMCTION  OF  z 
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Beglnnlng  with  the  sixth  x-plane,  the  averaging  was  eliminated  (the 
standard  algorithm  was  employed).  At  the  eleventh  x-plana  a geometric 
mesh-stretchlng  algorithm  of  the  following  form  was  used; 

, n ■ 11,12,13,...  . 

No  such  stretching  has  baen  employed  In  the  z-dlr action,  but  In  the 
future  It  may  also  be  required  for  rapidly  changing  profiles.  It  should 
be  noted  that  our  averaging  algorithm  was  required  In  both  the  x and  z 
directions  to  remove  all  oscillations. 

A mesh  refinement  algorithm  is  used  which  adds  or  delstes  points 
depending  on  the  relative  local  variation  In  the  truncation  error  of  the 
difference  equations.  Roughly  80  grid  points  in  the  n-dlrection  and 
11  grid  points  In  the  z-dlrectlon  are  employed. 

5.0  Results 

Computations  ware  performed  on  the  Berkeley  CDC  7600  machine.  A 
typical  calculation  required  about  6 minutes  of  CPU  time.  Figure  5 
indicates  the  external  and  Initial  velocity  distributions  which  have  been 
used  aa  a baals  for  our  calculations.  Tha  parameter  6 was  Introduced  as 
shown  to  vary  the  Initial  cross  flow  while  keeping  the  total  velocity 
constant  as  a rough  simulation  of  a fixed  supply  of  engine  mass  flow, 
the  velocity  profile  was  selected  to  have  a characteristic  fully  developed 
character  associated  with  turbulent  wall  Jet  flows.  Also,  the  surface 
shape  la  chosen  planer  for  simplicity  In  ths  ensuing  analysis.  Future 
aspects  of  this  effort  should  consider  the  "eating  up"  of  the  potential 
core  which  la  asaumed  to  occur  upstream  of  the  Initial  station  of  thla 
analysis.  The  parameters  and  were  choeen  to  provide  alupa  end 
value  continuity  of  the  profile  at  y ■ y . For  y > y , the  profile 

TRAX  niaX 

has  a half  Gauaslan  character  aasoclated  with  a free  jet.  For  y < y 

tnAX 

thft  profile  has  a boundary  layer  character#  In  tha  examples,  the  u and 
w initial  profiles  were  assumed  to  be  identical.  Moreover,  the  8 


i 
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distribution  wsa  salectad  to  be  qualitatively  similar  to  that  observed 
by  rake  surveys  representative  of  the  XFV-12A.  The  zero  cross  flow  case 
was  achieved  by  setting  C to  10 

Figure  6 dsmonstratea  decay  of  the  peak  velocity  with  the  standardised 
distributions  of  Fig.  5,  with  and  without  cross  flow.  It  is  evident  from 
the  figure  that  initial  cross  flow  has  a drataatlc  affect  on  enhancing  the 
decay  of  the  tuxlnum  velocity.  In  the  calculations,  tha  decay  exponent  n 
in  the  axternal  velocities  Is  aasuned  aa  1/2,  roughly  In  accord  with  a 
value  obtained  from  a two-dlaanslonel  line  alnk  aimulatlng  Inflow 
originally  propoaad  by  G.I.  Taylor. Both  atraamwlss  Incraaalng  and 
decreasing  croas  flow  caiss  are  shown.  In  tha  11ns  sink  modal,  ths 
Inflow  at  each  poaition  x ■ Xq  along  ths  jat  boundary  Is  datsmlnad  by 
strsaowlss  rats  of  change  of  mass  flux  according  to  conisrvatlon  of  mass 
applied  on  a rectangular  control  surfaca  in  tha  Jat  layer.  From  an 
observation  point  P In  the  axtamal  flow  which  la  aaaumad  Invlacld  and 
irrotational,  the  velocity  potential  can  therefore  be  rapraaentad  as 
an  isolatad  alnk  whoaa  intensity  is  proportional  to  ths  Inflow  at 
X * Xq.  Tha  cumulaclvs  affect  at  P of  all  such  sinks  at  x ■ x^  is 
obtained  by  a auparpoaltion  integral  of  all  thaaa  contributions  giving 
a llna  sink  raprasantatlon.  In  a more  raaliatic  modal,  thasa  axtamal 
vslocity  dlatrlbutlons  should  be  corrected  for  three>dlnensionallty  and 
alllptlc  Interaction  with  the  wall  Jet.  k calculation  of  this  typa  would 
be  a more  accurate  rapraaentatlon  than  the  present  approach  of  planforra 
and  surface  curvature  af facta.  Thaae  developments  are  strongly  recommended 
for  future  application.  In  this  connection,  we  recognize  that 
the  present  means  of  simulating  taper,  awsepback,  and  apanwlae  preaeute 
gradients  is  solaly  through  cross  flow  adjustment. 

Tha  thraa-dlmensional  invlacld  potential  ip  can  be  chsractarlzed  by 
surface  sink  distribution  of  ths  form  (see  Fig.  7) 


(x-t)*  -f  y*  + (z-C)* 
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where  S Che  area  otf  Integration  refers  to  the  total  Jet  area  on  and  off  the  wing. 
The  quantity  0 ia  the  sink  strength  obtained  by  matching  with  an  "outer  limit" 
of  the  second  order  solution  for  the  velocity  normal  to  the  body  appearing  in 
the  viscous  inner  wall  Jet  solution.  Generalizing  the  previous  concepts,  Eq.  (29) 
la  obtained  by  ouperposltion  of  a surface  distribution  of  elementary  sink 
lolutlona  of  the  form  (>  • where  R*  - (x-d)*  + y*  +•  (*-!;)*•  Note  that 
the  surface  distributions  do  not  Interact  in  the  sense  that  they  produce 
no  at  locations  other  than  their  own  Hence  (|)^  'w  local  Inflow 

analogous  to  that  praviously  dascrlbed  for  the  two-dimensional  case. 

The  quantity  o for  two-dimensional  boundary  layers  la  analogous  to  the 
straamwlaa  gradient  of  the  displacement  thickness  6'(x).  To  include 
lifting  surface  effects,  a surface  doublet  or  vortex  distribution  should 
be  added  to  Eq.  (29).  The  local  vortex  strength  can  also  be  determined 
by  matching. 

The  Inflow  velocity  related  to  the  sink  Intensity  a in  Eq.  (29)  is 
in  turn  a function  of  the  entreinment.  This  quantity  is  elao  significant 
from  tha  standpoint  of  the  tradeoff  between  skin  friction,  BLC,  and  rapid 
acceletatlon  of  the  secondary  in  compact  thraa-dimenalonal  thrust  augmantlng 
ejectors  such  as  thoss  employed  on  the  XFV-12A  and  upper  surfeca  bloving. 

In  Fig.  8,  tha  comparieon  between  cross  flow  end  the  absence  of  it 
gives  the  indicated  nominal  entrainment  veristlons  with  atraamwisa 
distance  where  nominal  entrainment  Q is  defined  in  the  figure.  In 
spite  of  the  oppreclable  increase  in  decay  of  the  maximum  value  of  u 
shown  in  Fig.  6,  and  cssultant  shaar  stress  in  Fig.  9,  only  a slight 
difference  in  entrainment  quantity  and  rate  is  shown  in  Fig.  8.  The 
difference  in  maximum  velocities  which  are  similar  for  v,  the  spanwlae 
component,  arc  presumably  related  to  the  enhanced  dissipation  associated 
with  cross  flow  and  that  implied  by  tha  eddy  viscosity  model.  The  lack 
of  a corrsaponding  increase  in  entrainment  rate  may  be  dua  to  nonlinear 
compensating  nf facta  built  into  the  turbulence  model  and  cannot  be  raadily 
explained  on  an  intuitive  beels  at  this  time.  In  this  connection,  other 
calculatlona  should  be  performed  for  which  the  streamwise  component  of  the 
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initial  velocity  is  held  fixed  rather  than  its  overall  magnitude  on 

introduction  of  cross  flow.  Furthermore,  trends  Involving  an  Increase 

in  decay  rate  of  peak  velocity  with  an  Increase  in  etitrainment  rate  that 

occur  in  two-dimensional  fret  jets  based  on  conservation  laws  and  similitude 

must  be  reassessed  for  non-eimllar  three-dimensional  wall  jeta  such  as 

those  considered  herein.  Here,  the  variable  etraamline  direction  through 

the  Jet  layer  must  be  considered  as  wall  as  the  peak  of  the  resultant 

velocity  rather  then  those  of  the  individual  velocity  components. 

If  u is  the  local  streamline  direction  on  the  wall,  under  certain 
3q 

circumstancee,  may  become  more  negative  with  increasing  Q.  It 
should  be  noted  that  the  axpression  for  entrainment  Q given  in  Fig.  8 


aaeumes  that  w 0 at  the  tip  2 - z, 


If  this  ia  not  the  case,  an 


additional  term  must  be  added  to  this  relation.  As  in  Figs.  9 and  10, 

1/2 

qualitatively  olmilar  behavior  ia  obtained  for  the  cega  in  which  w x . 

Associated  with  the  previous  results,  Fig,  10  shows  the  effect  of 
cross  flow  on  Jst  spreading  rats  relatsd  to  As  previously,  only 

smell  diffsrencea  are  indicated  for  the  cases  considered.  In  Fig.  11, 
however,  an  important  upstream  movement  of  tho  nominal  itparatlon  line  is 
Indicated  with  the  introduction  In  croas  flow.  Hare,  nominal  separation 

is  defined  to  occur  where  t ■ 0.  A more  pertinent  definition  li  ■ 0 

where  q is  the  veloultv  component  normal  to  an  envelope  of  the  surface 
n 

etreamllnes.  Implementation  of  the  latter  definition  is  nnvlaioned  in 
follow-on  studies  involving  primary/ secondary  entrainment  interactions. 
This  result  is  significant  with  respect  to  penalties  associated  with  taper 
and  sweep  in  three-dimensional  ejector  dlffusora. 

In  Fig.  12)  another  Important  consequence  of  cross  flow  is  examined 
in  connection  with  the  surface  streamline  pattern.  In  the  figure,  two 
casts  are  compared  involving  differing  amounts  of  cross  flow.  Lines  such 
as  AB  and  A'B'  represent  typical  streamlines  for  the  two  different  cross 
flow  cases  in  which  the  downstream  direction  Lb  in  the  direction  of  the 
arrows  In  the  rectangular  area  oorrusponding  roughly  to  the  planform  of  a 
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rectangular  wing.  In  this  interpretation,  OQ  would  be  the  leading  edge, 

RS  the  trailing  edge  and  QR  its  tip.  Significant  enhancement  in  tlownstream 

streamtube  contraction  is  obvious  with  increase  in  cross  flow.  This 

contraction  could  presumably  lead  to  end  wall  separation  of  the  type 

observed  in  typical  short  ejectors.  In  Fig.  13,  a similar  picture  Is 

1/2 

indicated  for  the  Increasing  w x cross  flow  case.  Although  some 
qualitative  similarity  exists  for  the  decaying  case,  it  la  evident  that 
the  possibility  of  a separation  surface  streamline  envelope  exists  In 
the  pattern  of  the  streamlines. 

6.0  Concluslona 

A class  of  cases  was  Investigated  roughly  possessing  initial  flow 
angularity  and  adverse  pressure  gradients  protntypic  of  those  on  the 
blown  surfaces  of  typical  propulsive  lift  systems  such  as  the  Navy/ 

Rockwell  XFV-12A  thrust  augmented  wing.  Results  obtained  from  the 
computational  model  Indicate  that  if  the  Initial  total  velocity  is  kept 
fixed,  then  the  introduction  of  the  cross  flow  enhances  the  freestream 
decay  rate  of  the  peak  of  the  velocity  component  in  the  freestream 
direction.  In  addition,  the  entrainment  quantity  and  its  rate  decrease 
with  increased  cross  flow.  The  three-dimensional  phenomena  not  only 
influence  the  effect  of  taper  on  the  boundary  layer  control  characteristics 
of  a Coanda  flap,  but  also  Indlc.nte  a "jet  shrink"  which  could  be  a 
mechanism  promoting  end-wall  separation.  To  our  knowledge,  our  model 
is  the  first  to  quantify  such  trends.  Both  should  be  considered  in 
the  design  of  any  propulsive  lift  system.  Finally,  the  effect  on  the 
prescribed  external  adverse  pressure  gradient  in  the  presence  and 
absence  of  cross  flow  has  also  been  examined.  From  the  limited  results, 
the  spanwise  separation  line  moves  progressively  further  upstream  with 
Increasing  cross  flow. 
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PART  11.  USERS  MANUAL 


Although  th«  computor  progrtm  givan  In  tha  Appendix  waa  daalgnad  to 
run  on  the  Lawranna  Barkalay  Laboratory  CDC  7600,  it  may  be  run  on  other 
CDC  6600  or  CDC  7600  computera  by  either  removing  the  plotting  routlnaa 
from  the  MAIN  program  or  by  adapting  tham  to  meat  the  requirements  of 
any  othar  facility. 


1.  Deck  Setup  (Berkeley  CDC  7600  with  plots) 

Job  Card 

anotapes 

RUN76. 

LGO,NL-60000. 

DISPOSE, PLOT-PL,R«[SEND  PLOTS  TO  W.D.  MURPHY/ SCIENCE 
•fCENTER/ROCKTOLl  INTERNATIONAL/ 

+P.O.  BOX  1083/ 

+THOUSAND  OAKS,  CA.  91360],  SC-BKY. 

7/8/9 

Source  program 
7/8/9 

CASE  7 UE-1.0/SqRT(X)  WE-1.0E-5'»Z/SQRT(X)  FUT  PLATE 
$ INPUTS  XSUPLY- . TRUE.  ,YSUPLY- . TRUE. , ZSUPLY-  . TRUE . .OPTPT-.  TRUE . S 

t 

second  column 
6/7/8/9 


Note  chat  the  DISPOSE  card  merely  tells  Berkeley  where  to  send  the  plots 
and  may  be  changed  for  othar  users.  Alao,  [ ■ 7/8  punch  and  ] ■ 0/2/8 
punch.  For  other  computer  facilities  the  RUN76  card  may  be  replaced  with 
FTN,  and  the  *NOTAPES  and  DISPOSE  cards  should  be  removed. 
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2-  Estimate  of  Running  Time 

Ths  aampla  cate  (tea  tha  next  acction  and  the  Appendix)  required 


389  aecomAa  on  the  Berkeley  GDC  7600.  In  general,  tha  running  vlll 
depend  upon  the  grid  aelected. 


CD  Berkeley  GDC  7600  (apaclal  feature  ••  Calcomp  plotter). 

Cil)  Any  other  GDC  6600,  GDC  7600,  CDG  175,  or  GDC  176  may  be  ueed 
If  tha  plotting  routlnaa  are  removed  from  the  MAIN  program. 


4,  Name  and  Level  of  Proarammina  Languaae  Deed  In  Proaram 
FORTRAN  IV. 


A.  Input-Output  Information 


Gloaaary  of  Input  Parametera 

XSUPLY  Logical  variable  which  la  .TRUE,  if  uanr  auppllee  the  atreeawiae 
meah  (x-meah) . Default  value  ■■  .FALSE.. 


YSUFLY  Logical  variable  which  la  .TRUE,  if  user  auppllee  tha  n-meeh. 
Default  value  ■ .FALSE.. 


ZSUPLY  Logical  variable  which  la  .TRUE.  If  user  auppllee  the  z-neeh. 
Default  value  - .FALSE.. 
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XA  Real  variable.  The  starcing  atreamwiae  station.  Default 

value  ■>  0.0.  XA  say  be  initialized  In  subroutine  XMGSH  if 
XSUPLY-.TRUE.. 

XB  Real  variable.  The  last  streanwiae  stacion.  Default  value  ■ 1. 

XB  may  be  inltlallieed  In  eubroutlna  XMESH  If  XSUPLY-.TRUE.  . 

YA  Real  variable.  The  left  end-point  of  the  n-niesh.  Default 

value  " 0.0. 

YB  Real  variable.  The  right  end-point  of  the  n-mesh.  Default 

value  " 20.  YB  nay  be  Inltiallaed  In  subroutine  YMBSM  If 
YSUPLY-.TRUE.. 

ZA  Real  variable.  The  left  end-point  of  the  a-neah.  Default 

value  ■ 0.0. 

ZB  Real  variable.  The  right  end-point  of  the  z-mesh.  Default 

value  « 1.  ZB  nay  be  Initialized  in  subroutine  ZKESH  If 
ZSUPLY-.TRUE.. 

lIXS  Real  variable  initializing  the  x-neah  In  setup  of  the  streamwlse 

mesh  and  ueed  only  If  XSUFLY- .FALSE. . Default  value  ■ l.E-S. 

FACX  Real  varlabls  representing  the  multiplication  factor  In  the 

setup  of  the  x-mesh  If  XSUPLY- , FALSE . , l.e. , X(N)-FACX*X(N-1) . 
for  the  first  few  points.  Default  value  ■ 1.2. 
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NX  Int«g«r  virlabU.  Th«  U*t  itriianiwii*  station.  This  may  be  sup- 

plied In  lubroutlne  XMESH  if  XSUPLY-.TRUE. . Default  value  - 60. 

J Integer  variable.  Total  number  of  ri-polnea  which  equala  the 

number  of  Internal  Intervale  plua  1.  J nay  be  supplied  by  uaer 
In  aubroutlne  YMESH  If  YSmV-.TSUE..  Default  value  - 101. 

I Integer  variable.  Total  number  of  z-polnta.  I may  bo  auppliod 

by  uaer  In  aubroutlne  ZMESH  If  ZSUPLY>.TRUE. . Default  value  ••  11. 

The  maximum  value  la  alao  11  unlaea  the  COMMON  BLOCKS  In  the 
MAIN  program  are  enlarged. 

REFINE  Logical  variable  which  la  .TRUE.  If  the  n-neih  la  to  be  refined. 

Default  value  - .TRUE. . 

4 

KMAX  Real  variable.  Maximum  n-meah  Interval  permitted  in  refinement 
routine.  Default  value  - 2.1. 

OPTPT  Logical  variable  which  la  true  If  uaer  auppllaa  atreamwlae 

atatlona  to  be  printed  on  paper.  These  will  alto  be  plotted 
if  the  Berkeley  7600  plot  routines  are  being  employed.  Default 
value  .FALSE..  Ueer  auppllee  these  stationa  in  subroutine 
PRMESH. 

1 

1 
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MPRINT  Integer  verlable.  Number  of  uniform  streamwlae  stutlone  at 
which  lolutlon  la  to  be  printed  ((XB-XA)/NPR1NT) . Default 
value  * 10. 

USUPLY  Logical  variable  which  la  .TRUE,  if  user  auppllaa  the  initial 

velocity  profllea  in  aubroutlna  PROFLE.  Default  valua  ■ .TRUE. . 

If  thaae  profllea  era  not  auppllad,  aimllarlty  aolutlona  will 
be  generated  and  uaad  aa  atarting  condltlona. 

R Real  variable.  Raynolda  number.  Default  valua  ■ 1.0/14.216862S-5 . 

Input 

Ttte  flrac  card  In  the  Input  atraam  containa  a TITLE  card  F0RMAT(8A10) 
daacrlblng  the  problem  to  be  aolvad.  The  remaining  carda  contain  NAMELIST 
data  of  the  varlablea  In  the  above  gloaaary.  The  Hat  la  defined  by 

NAMELIST  /INPUTS/  XSUPLY,YSUPLY,ZSUPLY,XA 

The  uaer  muat  alao  aupply  the  rouclnea  XMESH,  YMESH,  EMESH,  PRMESH, 
and  PROFLE  If  the  varlablea  XSUPLY,  YSUPLY,  ZSUPLY,  OPTPT,  and  USUPLY 
have  the  value  .TRUE.,  reapectively.  How  these  subroutines  are  written 
la  fully  deacrlbtd  in  the  program  Hating  In  the  Appendix. 

In  addition,  the  uaer  muat  aupply  the  routine  PREPP  which  defines  the 
values  of  P^,  P2,  ...,  P^^q  (equation  (14b))  In  terms  of  functions  of  u^, 
w^,  hj^,  h^,  K^,  and  K2.  In  the  example  given  in  the  Appendix,  •• 

■ 10”^  *//k,  hj^  ■ hj  ■ 1,  and 
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Plots 

Th«  progrsm  was  dsalgnad  to  run  on  tha  Barkalsy  CDC  7600  and  to  uas 
ths  Calcomp  plottar  subroutlna  packaga  avalltblt  at  Barkalay.  It  may  ba 
run  on  othar  CDC  eomputars  by  althar  ramovlng  tha  calls  to  CCNEXTi  CCGRIDi 
CCLTR,  CCFLOT  and  CCEND  in  tha  MAIN  program  or  by  adapting  tha  plot 
routinaa  to  maat  tha  raquirameuts  of  tha  naw  facility. 

Graphs  ara  aada  of  tha  following  quantltiaai 

1.  y varsua  u/u^  at  Sg,  and 

2.  y varsua  w/w^  at  a^,  a^,  Zg,  and  Zj^^^ 

3.  Shaar  atrass  varsua  x 

4.  Straam  linas 

5.  varsus  x 

6.  Jat  spraading  varsus  x 

7.  Entrainmant  quantity  varsus  x. 

Tha  labals  for  tha  individual  plots  ars  codad  in  LABELl,  LABEL2, 

....  LABEL7  of  tha  MAIN  program  and  may  bo  changad  for  each  individual 
case  study,  For  thasa  plots  tha  value  of  u must  be  codad  for  aach  naw 
case  in  the  MAIN  program  and  in  tha  subroutlna  NEVTON  whara  is  also 
raquirad.  In  tha  axampla  in  the  Appendix  ■ l//x  and  w^  ■ lO”^  a//x. 
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Example  (sample  input  deck)  ' 

CASE  7 UE-1.0/SQRT(X)  WE-1.0E-5*Z/SQRT(X)  FLAT  PUTE  ,,  ■ 

/ l' ' 

9INPUTS  X8UPLY- , TRUE . , YSUPLY- . TRUE . , ZSUPLY- . TRUE . , 

OPTPT-.TRUE.$ 

The  first  card  la  the  title  card  which  la  printed  on  paper,  and  the  aacond 
card  la  a NAMELIST  card  which  impllaa  that  chu  uaar  will  supply  his  own 
X'taash,  n-mosh,  a-maah,  and  the  straanwlae  stations  to  be  printed  on 
papar.  Examples  of  these  subroutines  can  be  found  in  the  Appendix.  All 
other  variables  will  be  assigned  their  default  values. 

The  initial  velocity  profiles  are  coded  in  subroutine  PROFLE  and 
represent  the  Initial  conditions  given  on  page  18  of  Part  I of  this 
report.  These  initial  profiles  may  easily  be  changed  by  following  the 
initructiona  given  in  the  subroutine. 


Output 

The  title  card  la  printed  at  the  top  of  the  paga.  If  REF INE^ .FALSE. , 
thn  grid  points  (x-ri“t)  and  y (assuming  u^-1)  era  printed.  This  is 
followed  by  the  Initial  profiles  where  F»f,  DF»f',  DDF^f",  G"g,  DG"g' , 
and  DGG-g".  The  indices  IZ  end  L represent  tha  subscripts  on  the  i- 
and  ri'mesh  pointa,  reapactivaly . If  REF INE". TRUE. , tha  program  refines 
the  Iniclel  mesh  and  prints  the  profiles  on  ths  naw  mesh  before  printing 


it  i ' 


the  refined  grid.  The  NAMELIST  le  printed  next. 


f 


/ V.  - . 
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Every  time  Newton's  method  converges  for  a given  streamwise  station 
and  2-polnt  the  following  daca  ia  recorded: 

1.  the  number  of  iterationi  Cor  convergence  (ITER) 

2.  maximum  absolute  error  (ERROR) 

3.  element  of  the  solution  vector  where  the  maximum  error  occur a 

4.  y*  defined  on  page  11  of  Part  1 (YTC) 

5.  f"(0)  (U(3,l)) 

6.  a"(0)  (T) 

^10'  '"20-  ’’30*  ^0*  So*  ^ 

8,  ahear  streaa  (TAU) . 

After  a oomplete  z-pUne  ia  swept  the  following  data  ia  printed) 

1,  entrainment  quantity  (Q),  page  24  of  Part  1 

“ ‘l-  ‘3'  ‘6'  ‘S'  ‘10 

“ ‘l'  ‘3'  's'  ‘S'  ‘lO' 

Note  that  u and  y are  defined  on  page  27  of  Fart  I. 
nidx  niAX 

The  complete  solution  is  printed  at  those  streamwlae  atatlona 
prevlouely  specified  by  the  user  in  PRMESH  or  at  default  values.  The 
format  for  thle  printing  is  the  seme  as  that  for  the  initial  profiles. 
Plot  data  will  also  be  loaded  into  storage  vectors  at  these  stations  for 
possible  later  use  by  the  plotting  subroutines  in  the  MAIN  program. 
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Example  (sample  output  listing) 

Baoauae  of  the  three-dimenalonal  nature  of  this  problem,  a aompleta 
listing  of  tha  output  would  be  too  axhaustlva  for  this  raport.  Therefore, 
we  give  below  the  first  page  of  output  that  raaulted  from  tha  sample  Input 
deok  given  earlier. 
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B.  Subroutine  Description 

All  subroutines  contain  many  comment  cards  describing  the  routine 
and  Its  function.  Our  purpose  here  Is  not  to  dupllrate  all  this  Informt- 

I 

tlon  but  Instead  to  discuss  some  of  the  technical  details  of  each  routine. 

An  * denotes  the  routl.ie  is  user  supplied. 

SUBROUTINE  NAME:  BC 

PURPOSE:  This  routine  computes  the  boundary  conditions  and  the  Jacobian 
of  the  boundary  data.  See  equation  (13)  of  Part  I. 

DESCRIPTION;  The  solution  vector  Is  vnrltten  In  the  form  u ■ (i^tf ' if" >gtg' >g")  i 
and  the  variables  A and  B denote  left  and  right  boundary,  respectively. 

The  boundary  conditions  are  assumed  to  be  written  In  homogeneous  form  and 
the  left  boundary  conditions  are  coded  before  the  right  ones.  Por  example, 
the  right  boundary  condition 

f'(nj  - 1 j 

translates  into  code  as  | 


i 


i 


I' 


G(5)  - UB(2)-1.0. 

SUBROUTINE  NAME:  BETASV 

PURPOSE:  This  subroutine  solvee  for  S in  the  LU-'decomposltlon  of  the  block 
trldlagonal  Jacob  lan  matrix  (equation  (28a)  in  Part  I) . 
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J" 
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DESCRIPTIONt  Tha  block  trldlagonal  matrix  Is  written  as 


/A  S <1-1,2 J) 


where  A^,  and  are  matrices  of  order  6. 


/A  la  decomposed  as 


/A  2 [6^  I OltO  C^] 


where 


1 . “l  ■ *1 

(la) 

1 

1 Sj  - Bj 

li'*2  fBye  te^hJ 

(lb) 

1 

i “l  ■ *1  - ®1  'l-l 

l.“2  ^SfeseyvT  • 

(Ic) 

For  simplicity!  we  write  the  LU  decomposition  of  as 


“l  ’ ^ “l  • 


(2) 


If  this  decomposition  doesn't  exist,  then  permutation  matrices  for 
row  and  column  pivoting  can  be  Introduced  In  equation  (2)i  Assuming  p 
rows  In  B^,  the  solution  of  In  equation  (lb)  Is  given  by 


T T 

’*1-1  ?k  - ?1 


T aT 
1-1  ?k  " ^k 


k-l,2,...,p  . 


(3) 
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SUBROUTINE  NAME:  BLOCK! 

PURPOSE i Thli  subroutine  decomposes  the  block  trldlagonal  Jacobian  matrix 
Into  LU-£ormi 

DESCRIPTION:  The  steps  in  equations  (1-3)  above  era  carried  out  by 
performing  the  indloated  matrix  multiplications  and  calling  routines 
LUSOLV  and  BETASV. 

SUBROUTINE  NAME:  BL0CK2 

PURPOSE:  This  routine  solves  the  block  trldlagonal  syatsm  /Ax  ■ b 
assuming  the  factorised  form, 

DESCRIPTION:  Write  z = (Zj^,  z^,  Zj)^  and  x 2 (Xj^,  Xj.  Xj)'^ 

T T 

where  ■■  and  Xj^  ■■  ^21*  •***  * Then 

BLOCK2  first  solves 

“ bi  (^*) 

5k  “ ~k  " 5k-l  k-2,3,...,J  (^b) 

and  then 

°‘j  *J  “ !j  ^5a) 

oij^  ^ ^ 1 "*  1 p g 1 1 s 1 2 • (51d ) 


Equation  (5a)  la  solvad  by  calling  USOLVE. 
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SUBROUTINE  NAMEt  BOX 

PURPOSE:  Thli  lubrouclne  lati  up  tha  block  crldlagonal  ayacam  of  aquatlona 
for  tha  attachnanc  llna  aquation  (lS-16)  of  Fart  I whan  KASE*i2  and  tha  3-D 
wall  jat  aquatlona  (23-27)  of  Part  1 whan  KASE«3. 

DESCRIPTION:  Tha  alananta  of  [B^,A^»C^]  ara  loadad  by  calling  tha  routlnaa 
BCt  RHSF)  and  RHSF2D  which  contain  Jacobian  Information  for  the  boundary 
condltlona  and  tha  varloua  dlffaranca  aquations.  Pravloua  station  infonua- 
tion  and  updating  of  tha  turbulanca  modal  la  inc.ludad  by  Incorporating  tha 
G matrix  (aia  Appandlx  for  tha  componant  daacrlptlon) . Finally,  tha 
aquatlona  ara  raarrangad  to  anaura  tha  fir at  diagonal  block  la  nonalngular. 

SUBROUTINE  NAME:  LUSOLV 

PURPOSE:  Thia  aubroutina  daccmpoaaa  a scalar  matrix  Into  LU-form  ualng 
a mixed  pivoting  atrategy. 

DESCRIPTION:  Sae  Analvaia  of  Numarical  Mathoda  by  E.  Isaacson  and 
H.B.  Kallar  for  a standard  diacuasion  of  Gausalan  alimlnation  (LU  dacom- 
positlon)  with  pivoting. 

SUBROUTINE  NAME:  MAIN 

PURPOSE:  Thia  la  tha  main  driving  program. 
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DESCRIPTION:  lb's  KV/o  program  Initializes  solution  vectors,  sets  default 
values,  reads  and  xfrlCas  NAMELIST  data,  and  calls  subroutines  that 
initialise  the  grid  and  the  velocity  profiles.  The  nuarchlng  schema  la 
oallad  from  MAIN  as  vail  as  all  plotting  packages. 

SUBROUTINE  NAME:  NETRUN 

purposes  This  subroutine  Is  used  to  add  or  delete  ri-meah  points  in  order 
to  improve  the  smoothness  of  the  solution  especially  in  regions  where 
large  boundary  layers  exist. 

DESCRIPTION:  The  local  truncation  error  of  the  numerical  scheme  Is 
approximated  at  the  mid-point  of  each  sub-interval  of  the  initial  grid. 

Points  are  added  or  deleted  so  that  this  error  will  remain  constant  on 
the  whole  interval.  NETRUN  is  only  called  for  the  first  x-polnt  and 
e-point,  If  the  user  requests  this  option  (REPINE* . TRUE .) . This  new 
improved  grid  ia  used  for  the  entire  3-D  region. 

SUBROUTINE  NAME:  NEWTON 

PURPOSE:  This  subroutine  solves  the  2-D  attachment  line  aquation  whan 
KASE*2  and  the  3-D  wall  Jet  equation  whan  KASE*3  by  employing  the 
undarrelaxad  Newton's  method  (equation  (28)  of  Part  I). 

I 

I 

i 
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DESCRIPTION:  Th*i  following  notation  la  uaad: 

N»N1  ■ numbar  of  unknown!  In  the  solution  vector  ■ 6 

NP"N2  ■>  number  of  boundary  conditions  deflnsd  at  the  left  boundary  ■ 4 

NQ"N3  * number  of  boundary  conditions  deflnsd  at  tha  right  boundary  2. 

Pravioua  station  information  is  loaded  from  subroutine  PREPG  and 
the  eddy-vlsooalty  from  PREPB.  After  the  proper  block  trldlagonal  Jacobian 
matrix  la  loaded  <KASE-2  or  3),  it  la  solved  uaing  BLOCKl  and  BL0CK2. 
Undarrelaxation  (in-. 7)  is  performed  on  this  solution,  and  a maximum 
abaoluts  or  relative  (REL*.TRUE.)  error  between  iterations  la  computed. 

If  this  error  la  greater  than  £FSERR»10~^  (loaded  from  MAIN),  the  above 
steps  are  repeated;  otherwise,  plot  and  print  vectors  are  loaded,  u^  and 
are  coded  into  this  subroutine  and  must  be  changed  for  each  problem. 

SUBROUTINE  NAME:  OUTPT 

PURPOSE:  This  subroutine  writes  the  solution  on  paper  for  Che  plane  X~XN. 

DESCRIPTION:  The  logic  Is  set  up  so  Chat  Che  initial  profile  will  only 
be  printed  once.  The  x-statlon  is  loaded  Into  tha  plot  vector  for  later 
use.  The  array  ULAST  contains  f and  g'  for  future  plots.  Tha  solution 
vector  Is  printed  on  paper  using  the  format  described  in  tha  Output  section. 
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SUBROUTINE  NAME:  PREP 

PURPOSE:  This  coutln*  up  thm  Initial  grid  and  loads  Initial 
condltlona  Into  tha  solution  vactor.  Tha  grid  la  prlntad  from  this 
subrout Inai 

DESCRIPTION:  Maah  points  usa  dafault  valuas  If  XSUPLY,  YSUPLY,  and 
ZSUPLY  ara  .FALSE..  For  aaoh  .TRUE,  valua  a usar  auppllsd  subroutlna 
Is  callad.  Valuas  of  x,  n«  a,  and  y ara  prlntad  whara  y la  eomputad 
froB  u^  at  som  x atatlon»  usually  X|i.  Tha  usar  should  chanaa  u^  with 
avary  naw  casu  study.  If  U8UPLY".TRUE. . Initial  conditions  ara  loadad 
fron  tha  uaar  suppllad  routina  PROFIE;  otharwlsa^  a similarity  solution 
la  ganaratad  fron  Initial  axponautlal  dacaylng  functions  glvan  by  stata- 
nant  functions  daflnad  by  F1(T),  F2(T)t  and  F3(T).  Grid  raflnanant  on 
althar  of  thaaa  initial  condltlona  will  ba  parfomad  If  REF1NE«.TRUE. . 

SUBROUTINE  NAME:  PREPB 

PURPOSE:  This  subroutlna  modala  turbulanoa  using  tha  two  layara  dlacusaad 
on  pagas  B through  11  of  Part  I. 

DESCRIPTION:  Tha  valuas  of  and  C2  (TCI  and  TC2,  raapactlvaly)  tiava 
baan  sat  to  saro  using  a DATA  STATEMENT  but  thay  nay  ba  asnlgnad  valuas 
batwaan  aaro  and  thraa.  Tha  progran  Is  dlvldad  Into  two  parts  dapandlng 
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upon  vhathar  ch*  attachmanC  line  equations  (IZ«1)  or  the  3-D  wall  jet 
equation*  (IZ***!)  are  being  aolvad.  In  aanh  part  tha  value  of  y*  (YTC) 
la  detamined)  and  tha  quantltiaa  b ■ In  equations  (11)  and  (12)  of 
Part  I are  conputad  and  atored  in  the  second  and  sixth  rows  of  tha 
a natrix,  raapactlvely.  Previoue  station  data  is  also  updated  and  atored 
in  tha  third  and  fourth  rows. 


SUBROUTINE  NAMEi  PREPG 


PURPOSE!  This  subroutine  computes  previous  station  data  for  tha  full 
3~D  vail  jet  equations. 


OSSCRZPTIONi  The  actual  variables  coded  in  this  routine  era  T (0(3, L)) 
and  S (G(4|L))  which  vara  first  derived  in  Ref.  9 of  Part  I.  The  following 
notation  is  used  for  the  FORTRAN  variables) 


uHx 


UHXX  - u' 


n-l,i 


uum  « u 


j-l/2 

n-l,i-l 


'j-1/2 
- “j-l/2  • 

SUBROUTINE  NAME)  PREPOZ 


PURPOSE)  This  subroutine  computes  previous  station  data  for  the  atCacIt- 
mant  line  aquation  (aquations  (13-16)  of  Part  I). 
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DESCRIPTION!  Th«  quantltla*  T^Ii/2  ®J-i/2 

Rtf.  9 tr«  oodtd  In  thlt  ■ubrouclnt.  UH  dtnoCtt  eh«  varitblt  Uj"}... 

SUBROUTINE  NAMEi  PREPP* 

PURPOSE!  Tha  vtrltbltt  P^,  P^.  ....  P^q  llvtn  by  aqumtlon  (14)  of 
Part  I tv*  oodtd  In  thit  tubroutlnt. 

DESCRIPTION:  For  aaoh  ntw  probltm  tha  followins  quant Itiat  mutt  ba 
codad  into  thla  aubroutlnat 

UK  ■ 

UEX  - 3u^/3x 
UEZ  • 3u^/3* 

WE  ■ 

WEX  ■ 3w^/3x 
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givan  by  aquation  (14)  • Nota  that  whan  tha  routlna  computaa  cha 

corraapondlng  valuta  of  P for  tha  attachnant  llna  aquatlona. 

SUBROUTINE  NAMEi  PRMESH* 

PURPOSE)  Thla  routlna  allowa  tha  usar  to  aupply  tha  atraanwlaa  atatlona 
at  which  aolutiona  art  to  ba  prlntad  on  papar.  Thaaa  atatlona  will  also 
larva  aa  thoaa  whara  plot  vaotora  will  ba  loaded. 

DESCRIPTION)  Saa  program  Hating  in  tha  Appandlx. 

SUBROUTINE  NAME)  PROFLE* 

PURPOSE)  Thla  routlna  aata  up  tha  Initial  valnolty  profllaa. 

DESCRIPTION)  In  tha  aaapla  caaa  glvan  in  tha  Appnndlx  wa  hava  codad  tha 
Initial  condition!  daacrlbad  on  paga  18  of  Part  I.  Tha  quantity  u /u 

BlAX  m 

varlaa  from  1.5  to  2.0.  Tha  flrat  alx  rowa  of  tha  USTORE  matrix  contain 

tha  valuaa  of  (f  .f  ,f"igig' .g")  aa  a function  of  rte  ‘nd  t . If 

^ in 

REFINE*. TRUE. , tha  maah  raflnaraant  routlna  altara  tha  n-maah  in  ordar 
to  maka  tha  Initial  conditions  amoothar  batwaan  grid  pointa.  Naw  Initial 
conditions  ara  cotnputad  and  prlntad  at  tha  rvfinad  pointa. 
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SUBROUTINE  NAME)  RHSF 

PURPOSE t Thin  toutin*  computas  tha  quantitlaa  on  tha  right-hand  a Ida  of 
tha  tquatlons  glvan  by  <31-36)  in  Raf.  9.  Tha  Jacobltn  nattlx  It  alao 
•valuattd  and  ttored  In  cha  A matrix. 

DESCRIPTION)  Tha  following  notation  la  uaad) 

UT  - u/4.0 

UBAR  ■ u 

U ■ currant  valua  of 
UX  ■ paai  valua  uj’ 

UXX  ■ paaa  valua 
UXXX  - paaa  valua 

Tha  valuta  of  tha  right-hand  alda  of  equations  (31-36)  are  storad  In 
F(l),  F(Z),  F(6),  raapectlvely 

SUBROUTINE  NAME)  RUSF2D 

PURPOSE)  RHSF2D  la  similar  to  RHSF  except  tha  attach))iant  Una  aquations 
glvan  by  (38)  In  Raf.  9 and  tha  corresponding  Jacobian  matrix  are  coded. 

DESCRIPTION)  Saa  RHSF  above. 
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SUBROUTINE  NAME!  TRUN 

PURPOSE!  Thla  lubroutin*  oomputta  th«  local  truncation  arror  of  tha 
oantarad-Eular  achama. 

DESCRIPTION!  Saa  H.B.  Kallar,  '*Aocurata  Dlffaranoa  Mathoda  for  Nonllnaar 
Two-point  Boundary  Valua  Problamat”  SIAM  J.  Numar.  Anal.,  11  (1974), 
pp.  303-320. 

SUBROUTINE  NAME!  USOLVE 

PURPOSE!  Thla  routlna  aolvaa  tha  scalar  matrix  aquation  Ax  ■ f aftar 
A has  baan  put  In  LU-fom. 

DESCRIPTION!  Saa  E.  Isaacson  and  H.B.  Kallar.  Analvals  of  Numarloal 
Mathods.  J>  Wllsy  & Sonst  Naw  Yorki  1966. 

SUBROUTINE  NAME!  VANDET 

PURPOSE!  Thla  function  oomputas  tha  datarult\ant  of  an  nxn  Vandannonda 
matrix  for  1 < n < 7.  This  routlna  Is  callad  by  TRUN  In  ordsr  to  computa 
tha  local  truncation  arror  of  tha  oantarsd-Eul > ; or  box  achama. 
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DESCRIPTION:  Thu  Vandermonda  dacannlnant  la  given  by 


n (x.  - xO  . 
i>J  ^ J 


SUBROUTINE  NAME!  VELHAX 

PURPOSE:  Thii  aubroutlna  computaa  and  y „ daicrlbad  on  page  18 

uUIX  iMlX 

of  Fart  I. 


DESCRIPTION:  vi(rij)  is  aomputad  for  each  meah  point  Hj  end  the  maxinum  | 

value  la  datarmlnedi  say  u ■ u(n*).  Thar  y le  given  by 

uuix  nuix 

y • r]*t/vx/\i  4 

'^max  ' a 

SUBROUTINE  NAME:  WALJET 

PURPOSE:  Thla  la  the  vain  marching  routine  which  calle  NEWTON  to  solve  I 

althar  the  3-D  well  Jet  problem  (aquatlona  (11-14)  of  Pert  I),  or  the  2-D 
attaohoant  llna  equations  (equations  (13-16)  of  Part  I). 
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DESCRIPTION!  Th«  marching  direction  1b  dlsouaBcd  on  page  16  of  Part  I. 
Unnatural  oBclllatloni  are  dampad  out  by  the  averaging  procaduro  diacuiead 
on  pagea  17  and  19  of  Part  I.  The  OUTPT  routine  la  callad  at  the  approprlata 
point  from  thla  aubroutlne  dapendlrg  upon  whether  OFTFT  la  iTRUE.  or 
.FALSE., 

If  aaparatlon  (f"(0)  < 0)  occura,  Che  marching  la  taminatad 
and  Che  aolutlon  at  the  pravloua  atreanwlaa  atatlon  la  printed. 

SUBROUTINE  NAME:  XMESH* 

PURPOSE!  Thla  routine  allowa  the  uaer  to  input  hia  own  x-maah. 

DESCRIPTION!  See  the  Appendix  for  an  axample  and  a conpleta  explanation. 
SUBROUTINE  NAME!  YMESH* 

PURPOSE!  Thla  routine  allowa  Che  uaer*  to  Input  hla  own  n-uash. 

DESCRIPTION:  See  the  Appendix  for  an  axample  and  a complete  explanation. 
SUBROUTINE  NAME:  ZMESH* 

PURPOSE:  Thla  routine  allowa  the  uaer  to  input  hla  own  z-maah. 

DESCRIPTION!  See  the  Appendix  for  an  example  and  a complete  explanation. 
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